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Highly frustrated antiferromagnets composed of magnetic rare-earth moments are currently at- 
tracting much experimental and theoretical interest. Rare-earth ions generally have small exchange 
interactions and large magnetic moments. This makes it necessary to understand in detail the role 
of long-range magnetic dipole-dipole interactions in these systems, in particular in the context of 
spin-spin correlations that develop in the paramagnetic phase, but are often unable to condense 
into a conventional long-range magnetic ordered phase. This scenario is most dramatically empha- 
sized in the frustrated pyrochlore antiferromagnet material Tb2Ti2 07 which does not order down 
to 50 mK despite an antiferromagnetic Curie- Weiss temperature Tew ~ — 20 K. In this paper we 
report results from mean-field theory calculations of the paramagnetic elastic neutron-scattering 
in highly frustrated magnetic systems with long-range dipole-dipole interactions, focusing on the 
Tb2Ti207 system. Modeling Tb2Ti207 as an antiferromagnetic (111) Ising pyrochlore, we find that 
the mean-field paramagnetic scattering is inconsistent with the experimentally observed results. 
Through simple symmetry arguments we demonstrate that the observed paramagnetic correlations 
in Tb2Ti2C>7 are precluded from being generated by any spin Hamiltonian that considers only Ising 
spins, but are qualitatively consistent with Heisenberg-like moments. Explicit calculations of the 
paramagnetic scattering pattern for both (111) Ising and Heisenberg models, which include finite 
single-ion anisotropy, support these claims. We offer suggestions for reconciling the need to restore 
spin isotropy with the Ising like structure suggested by the single-ion properties of Tb 3+ . 

PACS numbers: 75.10.Hk,75.25.+z,75.30.Gw,75.40.Cx 



I. INTRODUCTION 

The pyrochlore oxides, with the general formula 
(A2B2O7), have attracted a great deal of attention over 
the last decade because the combination of lattice ge- 
ometry and chemical composition allow for a plethora 
of interesting physical phenomena in these materialsi 2 ^. 
The A and B sites reside on two distinct interpenetrat- 
ing pyrochlore networks of corner sharing tetrahedra, 
Fig. ^ Since A 3+ can be either a magnetic rare-earth or 
a non-magnetic transition metal, and B 4+ can be a tran- 
sition metal with or without a moment, there are numer- 
ous possibilities to study insulating and itinerant mag- 
netic models in a geometrically frustrated environment. 
Experimentally, long-range magnetic ordered states^ 5 ^, 
novel magnetic phases (e.g., spin-glass£&2ii2iii, spin- 
ico 3 i 12 i 13 ' 14 , spin liquid 15 : 16 ), anomalous Hall effect, 
metallic properties^, and superconductivitjii2i2£ have 
been observed in the pyrochlore oxides. Materials with 
the pyrochlore-related spinel structure have also at- 
tracted much attention recently. Heavy fermion physics 
has been observed in the d-electron L1V2O4 compound 21 . 
In the spinel antiferromagnet ZnCr2F4, a spin-Peierls-likc 
transition has been observed^ at low temperatures and a 
protectorate of weakly interacting spin directors slightly 
above the spin-Peierls transition temperature 2 ^. 

In insulating rare-earth magnetic pyrochlores, the 
magnetic rare-earth ions often have large dipole mo- 
ments, i.e., (J, » 1/ib and small Heisenberg exchange 
interactions. In such a situation, long-range dipole-dipole 




FIG. 1: The pyrochlore lattice structure. The lattice is con- 
structed from a network of corner sharing tetrahedra in which 
each tetrahedron resides at an fee Bravais lattice point. A cu- 
bic cell contains 4 tetrahedra or 16 atoms. 



interactions are a significant contribution to the Hamil- 
tonian. As well, a large single-ion anisotropy interaction 
is also often present. For example, crystal fields pro- 
duce an effective Ising doublet at the rare-earth ion site, 
A 3 +, in the HoaTioOy (Ref. El, Dy2Ti 2 7 (R e f. |24|), 
Tb 2 Ti 2 7 (Refs.Hjandii an d Yb 2 Ti 2 7 (Ref.lH) ma- 
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terials. An energy gap separates the ground state dou- 
blet from the lowest lying excited states with the Ising 
quantization axis coinciding with the local cubic (111) 
directions'^, i.e., the quantization axis points toward 
the center of the tetrahedral basis unit cell. Most no- 
ticeably, this strong local (111) Ising axis anisotropy has 
been found responsible for endowing frustration to the 
pyrochlore lattice in the presence of effective nearest- 
neighbor ferromagnetic (FM) interactionsi2i2L2&. From 
a statistical mechanics point of view, a FM (111) Ising 
pyrochlore model is equivalent to a model for disordered 
water ice, It, where both magnetic^ and water iceSi^Sii^i 
models possess macroscopic degeneracy. Recent experi- 
ments on the ma gnet ic systems Dy2Ti207 (Ref. IT^ l and 
H02T12O7 fR.efs.l32l and IssT) . so called spin- ice mate- 
rials, reveal a residual entropy in agreement with the 
prediction for water-ice 31 . In contrast, nearest-neighbor 
antifcrromagnetic (AFM) interactions in a (111) Ising 
pyrochlore model are non-frustrated2£i2L2£i2i; therefore, 
such a model is expected to order at a temperature set by 
the nearest-neighbor energy scale, i.e., T N w J^/k B . A 
candidate AFM (111) pyrochlore material is Tb2Ti207, 
where a non-collinear long-range ordered state is ex- 
pected at Tn ~ IK (Refs. EH and l35|) . However, ex- 
periments indicate the material fails to order down to 
T = 50 mK (Refs. d and El), despite an antiferromag- 
nctic Curie- Weis temperature, 6*cw ~ —20 K, making 
this system so far one of the cleanest realizations of a 
spin liquid in a three dimensional system. The mecha- 
nism responsible for Tb2Ti207 failing to order down to 
such low temperature has not yet been resolved^. 

In Tb2Ti 2 07, the energy gap between the ground state 
doublet and the first excited state doublet is A w 20 
K (Refs. and 124^ . The size of anisotropy gap here 
is an order of magnitude smaller than what is observed 
in the spin- ice materials, A = 250 - 350 K (Ref. H3), 
and is comparable to #cw f° r Tb2Ti 2 07. Consequently, 
in early theoretical work on Tb2Ti2C>7 a (111) Ising 
model was assumed and a non-collinear long-range or- 
dered state with zero net moment about each unit tetra- 
hedron (q = state) was predicted at T « 1 K (Refs. I35I 
and l25|) . However, muon spin relaxation measurements 
indicate a dynamically fluctuating state at all experi- 
mentally achievable temperatures, which recent results 
pushed down below 50 mK (Ref. 16). In contrast, some 
static susceptibility data suggest a spin-glass state at 
these temperatures^. Despite the experimental evidence 
that the magnetic single-ion ground state of Tb 3+ is an 
Ising doubleliSiSi, several experimental groups have used 
simple Heisenberg type models (which for the pyrochlore 
lattice is paramagnetic down to T — > + , Refs. 38 39 40) 
to obtain qualitative agreement with the observed para- 
magnetic (PM) scattering patternii^. Finally, we note 
that high pressure neutron-scattering experiments on 
Tb2Ti2 07 find a transition to long-range order at a tem- 
perature in excess of 1 K for an applied pressure larger 
than 2 GPa^ 3 ., but the magnetic structure has not been 
determined. In summary, there is currently ambiguity 



even as to the nature of the paramagnetic correlations de- 
veloping in Tb2Ti207 at temperatures above #cw- This 
is an important issue since very recent neutron-scattering 
measurements on a single crystal of Tb2Ti207 down to 
50 mK have found that the scattered intensity in recipro- 
cal space remains essentially unchanged (i.e. frozen out) 
when going from the paramagnetic temperature of 10 K 
down to 50 mK fRef. Il6l) . It would therefore appear that 
a first requirement to make theoretical progress in un- 
derstanding the spin liquid state in Tb2Ti207 would be 
to understand the nature of the paramagnetic spin— spin 
correlations. It is the aim of this paper to shed some light 
on the paramagnetic correlations of this material. 

In this article, we use mean-field theory (MFT) to 
investigate the paramagnetic spin-spin correlations of 
a (111) Ising dipolar model pertinent to Tb2Ti207 as 
well as a finite Ising anisotropy model. In MFT, the 
PM regime is realized for temperatures above an en- 
ergy scale set by the mean-field (MF) critical mode, i.e., 
T > T C MF = A^ F /n, where T is temperature in units 
of fee, n is the number of spin components, and A^ F 
represents the mean-field global maximum eigenvalue of 
the q-dependent susceptibility. As the temperature ap- 
proaches T C MF , one expects the critical mode to control 
the spin-spin correlations. A clear understanding of this 
critical mode softening relies on an accurate treatment 
of all interactions in the Hamiltonian. The long-range 
dipole-dipole interactions are our main concern. Our ma- 
jor results are the following: We establish, on symmetry 
grounds and through calculations, that the observed PM 
neutron-scattering is inconsistent with a local (111) Ising 
dipolar model. Calculations performed for an anisotropic 
Heisenberg pyrochlore model yield good agreement with 
experiment and thus support the claim that at least a 
partial restoration of spin isotropy occurs in Tb2Ti207. 

We also have in mind a broader perspective in pre- 
senting the enclosed work and detailed derivation of 
the mean-field formulation of the structure factor, 5(q), 
for highly frustrated magnets. In the past few years, 
there have been a number of interesting and puzzling 
thermodynamic data and neutron-scattering results on 
highly frustrated magnets. Examples include the anti- 
ferromagnetic GdqGasO ^ garnet (referred to as GGG), 
(Refs. l44U45U46U47ll48M49ft the p yrochlore antiferromagnet 
Gd 2 Ti 2 07, (Refs. Hl5l l6l5Cl5lh and the pyrochlore ferro- 
magnet Yb2Ti207 (Refs. l2aand I52T) . The experimental 
results on these systems, which will be discussed in more 
detail in Section I1VI raise a common issue: To further 
our current understanding of a number of highly frus- 
trated magnets, including quantum fluctuations^, we 
need to have a clear and quantitative understanding of 
the predominant correlations that initially develop out of 
the paramagnetic state as the materials are cooled. The 
mean-field theory described herein is a first approach to 
formulate such a program. For concreteness, this paper 
focuses on the specific cases of insulating pyrochlore ox- 
ides with local (111) Ising axis anisotropy. It is straight- 
forward to use the formalism herein to tackle the frus- 
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trated dipolar systems described above as well as others. 

The outline of the paper is as follows: In Section ITU a 
condensed description of the MFT formalism of neutron- 
scattering for the anisotropic Heisenberg and (111) Ising 
pyrochlores is presented. We give the MFT results for 
the paramagnetic scattering of Tb2Ti2C>7 in Section ITTT1 
Our mean-field data are also compared to Monte Carlo 
results for the paramagnetic S(q). Section llVl discusses 
the need to relax the (111) Ising constraint to allow for 
transverse fluctuations in Tb 2 Ti 2 07. For completeness, 
we include several detailed appendices. In Appendix 1X1 
the variational MFT is discussed and the equations for 
elastic neutron-scattering are derived in detail. An alter- 
nate approach, a high-temperature series expansion, to 
the equations for neutron-scattering is presented in Ap- 
pendix [H] The derivation of the Ewald equations for the 
q-dependent dipole-dipole Hamiltonian is given in Ap- 
pendix |0 Appendix [D] contains a detailed discussion 
of the symmetry allowed scattering patterns for (111) 
Ising and Heisenberg pyrochlores. We note that other 
authors have discussed various aspects of the mean-field 
theor^r??!? 4 !^!^!^ and EwaldSMMMi derivations pre- 
sented in Appendices lAlandlUl respectively. Our purpose 
here is to provide a self-contained reference for the appli- 
cation of the Ewald method within a mean-field formal- 
ism, connecting at times with results from earlier work^, 
and which will be useful to other researchers who wish 
to study frustrated magnetic systems with non-Bravais 
lattice geometries and allowing for an array of possible 
spin symmetries and interactions. 



A. Models 

The pyrochlore lattice, Fig. [I] is a non-Bravais lat- 
tice that we describe as a fee lattice with a four atom 
unit cell. The positions of the fee Bravais lattice points, 
which coincide with a corner point on the tetrahedral 
basis, are denoted by Ri. The four atoms that form 
the tetrahedron at each fee point (and represent differ- 
ent sublattices) are labeled by r a . Hence, the position of 
a site in the pyrochlore lattice is given by Rf = Ri+r a . 
Table[I]lists our convention for the tetrahedral basis coor- 
dinates, r a . The most general Hamiltonian for rare-earth 
spins on the pyrochlore lattice is a Heisenberg model with 
nearest-neighbor exchange, dipole-dipole, and single-ion 
anisotropy (with a local (111) orientation) energies, 



H H = -J 



£ Sf ■ S b - ■ Sf) 2 

((i,a),(j,b)) 
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The unit vector z a represents the local (111) quantiza- 
tion axis that points toward the center of a tetrahedron. 
Table |U defines our convention for z a . The spins Sf have 
unit length and full 0(3) symmetry, Rf$ = Rf - R) is 
the vector separation between spins Sf and S b , and J 
and A define the exchange and single- ion energy scales, 
respectively. The convention established in Eq.^ defines 
J > as FM and J < as AFM exchange energies. The 
dipolar energy scale is set by Ddd = DR^ D , where 



II. MEAN-FIELD THEORY OF NEUTRON 
SCATTERING IN THE PYROCHLORES 



In this section, we present the models and provide 
an outline of the derivation of the mean-field equations 
for the neutron-scattering cross section, da(Q)/dfl, for 
classical spins on the pyrochlore lattice. The resulting 
equations are only applicable in the disordered param- 
agnetic regime of the model Hamiltonians. Our deriva- 
tion is performed for the general anisotropic Heisenberg 
model because it has broad appeal to the study of many 
highly frustrated magnetic systems. Our method is best 
described as variational mean-field theory^ (VMFT) 
(which is equivalent to a Gaussian approximation of the 
free energy) and has been used to study frustrated mag- 
netic systema 39 i 56 i 57 . The details of VMFT and its ap- 
plication to the scattering cross section are presented in 
Appendix [X] We remind the reader that MFT corre- 
sponds to a partial resummation of an infinite number 
of terms in the high-temperature series expansion for 
the q— dependent susceptibility, x(q). In particular, the 
correlations (S(0) ■ S(r)} are correctly treated to order 
(3 = 1/T (or x(q) to 1/T 2 ) in MFT. This is demonstrated 
in Appendix iBl where scattering cross section equations 
are derived via a high temperature series expansion. 
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H is the moment on the rare-earth ion, and R nn is the 
nearest neighbor distance. 

The (111) Ising dipolar model for the pyrochlore 
lattice^ is obtained by considering the limit of large Ising 
anisotropy in Eq.[TJ(A/| J| > 1 and A/D > 1). The low 
energy physics of this system is modeled by the Hamil- 
tonian, 



H, 



-J 
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By low energy physics we mean that the single-ion term 
in Eq. ^ is removed and the spins are restricted to lie 
along the local (111) quantization axis, i.e., Sf = z a af 
with af ± 1. If one were to truncate the dipolar sum in 
Eq. |21 at nearest neighbor distances, then the following 
effective nearest-neighbor energy scale could be defined, 



7 off = 7 



A. 



(3) 



where J nn = J/Z and D nn = 5D/3. For FM effective 
nearest-neighbor exchange, > setting all dipolar 
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TABLE I: Our convention for vectors: The r a define the basis 
vectors and z a define the local (111) anisotropy axes for spins 
on the pyrochlore lattice. The size of the corresponding cubic 
cell is given by a and contains 16 atoms. Vectors n w represent 
the global Cartesian basis vectors. 
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interactions beyond nearest-neighbor to zero, one has the 
nearest-neighbor spin-ice model of Harris et al, Refs. Il2l 
and I2J If the nearest-neighbor interactions are 
AFM, then the model possesses a unique ordered state 
(q = 0, all-in all-out state) at temperatures on the order 
of |J^| (Refs. l25l27l28l3# . The transition between the 
spin-ice and q = phases occurs at J nn /D nn = —1.0. 
When dipole-dipole sum is extended to long-range dis- 
tances, the transition between the q = and the spin- 
ice states shifts to J n n/D nn = —0.908. Hence, long- 
range dipolar interactions favor the q — AFM phase 
slightly^. 



B. Mean- field theory 

We are interested in calculating the elastic neutron- 
scattering cross section for both Heisenberg and (111) 
Ising spins on the pyrochlore lattice at the mean-field 
level. Therefore, we use the general anisotropic Hamilto- 
nian, Hh, as the starting point for MFT and include a 
local, fictitious field term, \h\ = \hf\, (where at the end 
of the calculation hf — > 0), 



i,j a.b u,v 



where 



(4) 



(5) 



Ad i<j d a ' b (z a -n u )(z b 



-D dd 



ab\ 



In the notation of this general model, the spin vectors are 
represented by Sf = ft* 1 ) S-' 1 +# S-' 2 Sf 3 , where 
the unit vectors n u are the global Cartesian unit vectors, 
see Table U and S°"' u is the u-th. component of spin. The 
sum in Eq. 0|does not include terms with Rf b = 0. For 
(111) Ising spins, one begins with Hi, Eq.|21 and adds the 
field term — J2i afof ' z a )crf. The resultant interaction 
parameter, J ah {i,j), does not include the spin compo- 
nents. 



The general expression for the mean-field free energy 
is as follows, 



F p = Tr{pH H }+TTr{pkip} 
= (H H ) p + T(lnp) p , 



(6) 



where p is the many-body density matrix and Tr repre- 
sents a trace over the states of p. A mean-field form for 
F p is obtained by first expressing the many-body density 
matrix as a product of single-particle density matrices 
p({S?}) = Ui, a Pi( S i), followed by minimizing F p with 
respect to pf (the variational parameters) subject to the 
constraints Tr{pf} = 1 and Tr{pfSf } — mf, where mf 
is the local, vector order parameter. For Ising spins, mf 
has only one component, mf. Next, the resulting mean- 
field free energy is transformed to momentum space by 
applying the definitions, 



q 

— £^(<z)e l *-<, (8) 



where N co u is the number of fee Bravais lattice points. 
We note that the above convention for the Fourier trans- 
form, which employs the position of the spin, Rf, results 
in a real symmetric q-dependent interaction matrix ^T(q), 
12 x 12 for Heisenberg and 4 x 4 for Ising spins. An al- 
ternate convention for the Fourier transform uses Ri , the 
Bravais lattice points, instead of Rf and yields a complex 
J{q), refer to Appendix lAl for details. For a non-Bravais 
lattice, the interaction matrix is not fully diagonalized 
by a Fourier transform. Hence, to completely diagonal- 
ize J{q) one must transform the q-dependent variables, 
m°, to normal mode variables. In component form, the 
normal mode transformation is given by 



a—l fi—1 



(9) 



where the indices (a, p) label the normal modes (12 for 
Heisenberg spins), and are the amplitudes of the 

normal modes. In matrix form, U(q) is the unitary ma- 
trix that diagonalizes J{q) in the spin£g)sublattice space 
with eigenvalues X(q). Hence, U^ p (q) represents the 
(a, u) component of the (a, p) eigenvector at q with 
eigenvalue X p (q). Finally, the mean- field free energy to 
quadratic order in the normal modes reads, 



A"(g))IC'f 



^E 



-q ' 



(10) 

where T P {T) = F P {T)/N ccih h?" oc h% u /T, T is the 
temperature in units of fee, and n — 3 for Heisenberg 
spins. Note, in order to consider the Ising case, the 
indices u and p are dropped from Eq. 1101 and n = 1. 
We have also dropped a constant from the expression for 
J- p {T), refer to Appendix lAl 
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The neutron-scattering cross section for unpolarized 
neutrons in the dipole approximation is given by the gen- 
eral expression 62 - 63 , 



with 



da(Q) _ C[f(Q)f 



N,, 



i,j a,b 



where Q is the momentum transfer, Q = G + q, Gisa 
reciprocal lattice vector and q is a vector in the first zone, 
f(Q) is the magnetic form factor of the relevant scatter- 
ing ion, and C is a constant. The spin-spin correlation 
function only involves spin components perpendicular to 
Q (i.e., S? ± = Sf - (Sf ■ Q)Q/\Q\ 2 ) and can be written 
as, 

( s t± ■ s j±) = E(" u • * v - ^ u ■ ■ Q)) 



x (S«' u S b j > v ) 



= E [ 5U,V - i^r) (srs b i v ) (12) 



where Q = Q/\Q\. The correlation function (Sf' u Sj' v ) is 
expressed as a thermal average of the mean- field variables 
and then transformed to normal modes, 

(SrS b S v ) =EEE(C^)% 8 (9)^(«') (13) 

q,ql a,f3 H,v 

x e - iq R h- tq, R 1 . 

The correlation function of normal mode variables is cal- 
culated from derivatives of the mean-field partition func- 
tion, Z = Trie-P^r^}, where T P (T) is given by Ea.Hni 
with respect to . The result is 



" } (3-A°(g)/T) 



for Heisenberg spins and 



t>q<l>qf) 



(1 - A"(g)/T) 



(14) 



(15) 



for (111) Ising spins. 

Using Eqs.^lEl andll4lorll5lin Ea. II II and carrying 
out the sums, one obtains equations for the scattering 
cross section. In the case of Heisenberg spins, we have 

Nccn d n - c[fm ^(a-A-O/T)' (16) 

where 



iGr a 



(17) 



is a three-component vector. For (111) Ising spins, one 
has 



N, 



1 da{Q) 
dQ 



F c l{q) 



cell 



c [/W)] 2 Ett^ 



(1-A«(qr)/T) 



(18) 



G r a 



(19) 



where F^{q) is still a three-component vector and z\ = 
z a - (z a ■ Q)Q. 

Equations 1161 and 1181 are the main results of this sec- 
tion. They provide a mean-field description for the PM 
elastic neutron-scattering of Heisenberg and (111) Ising 
moments, respectively, on the pyrochlore lattice. The 
temperature that defines the paramagnetic regime is set 
by the maximum eigenvalue according to 



T>T, 



MF _ 



max q {A max (q)}/r 



(20) 



where A max (q) is the maximum eigenvalue at wave vector 
q, and max q selects the global maximum for all q. The 
maxq{A max (q)} occurs at the ordering wave vector q or d- 



III. NEUTRON SCATTERING OF TB 2 TI 2 7 

Starting with the zeroth order (low energy) (111) de- 
scription for Tb 2 Ti 2 C>7 we use J\ = —2.64 K and 
D = 0.48 K yielding J nn /D nn = -1.1 (Ref. HJ, which 
compares to J nn /D Dn = -0.22 for Ho 2 Ti 2 07 (Ref. HI 
and J„„/Am = -0.52 for Dy 2 Ti 2 7 (Ref.|3f|. There- 
fore, at a nearest-neighbor cutoff distance, Tb 2 Ti 2 07 
is an AFM (111) pyrochlore that is predicted to de- 
velop non-collinear AFM order, with ordering wave vec- 
tor q old = 0, at T w 1 K (Refs. |H and HI • We em- 
phasize that in the context of a (111) Ising model with 
Jnn/Dnn = — 1.1, Tb 2 Ti 2 C>7 is still predicted to be a 
long-ranged AFM when both antiferromagnetic nearest- 
neighbor exchange and long range dipole interactions are 
considered 3 ^. Hence, the antiferromagnetic exchange in a 
(111) model of Tb 2 Ti 2 C>7 is sufficiently strong to prevent 
the perturbations arising from long-range dipolar inter- 
actions from changing the ordered state of the model. 
The counter point to the above model predictions is 
that experimentally Tb 2 Ti 2 C>7 remains a collective para- 
magnetic down to very low temperatures, T ~ 50 mK 
(Ref. El). 

As an initial attempt to explain the physics of 
Tb 2 Ti 2 C>7, we investigate the PM correlations within 
MFT and compare to the experimental results for elas- 
tic neutron-scattering. Experimental data for elastic 
neutron-scattering in Tb 2 Ti 2 07 are shown in Fig. Ufa). 
The most intense region of scattering is centered around 
Q = 0, 0, 2 with reduced correlations extending toward 
Q = 2, 2, and a scattering minimum at Q = 0. From 
the pyrochlore lattice structure and the MF formalism, 
we know that the intensity at Q — 0, 0, 2 is controlled 
by the eigenmodes at Q = 0, 0, but modulated by the 
phase factor exp(zG • r a ), see Eqs. 1751 and ITDl This 
raises the question as to whether the maximum about 
Q = 0,0,2 could be interpreted as the precursor of a 
long-range ordered non-collinear AFM state. 
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We begin by considering (111) Ising spins on the 
pyrochlore lattice. The details upon which our ar- 
guments here are based are provided in Appendix [D] 
The neutron-scattering intensity profile is determined 
by F"(q), Eq. 1191 and contains information on the 
spin anisotropy via z a and the eigenvalues A"(<jr), and 
the symmetry of the lattice through the eigenvectors, 
U a ' a {q), and a phase factor, exp(iG-r a ). Hence, the 
nature and strength of the exchange and dipole-dipole 
interactions are arbitrary. From these basic symmetry 
components, we find that |F£(0,0,0)| 2 = |Fj*(0, 0, 2)| 2 , 
or that the scattering intensity about Q = 0, 0, and 
Q = 0, 0, 2 has the same numerical value, disregarding 
the form factor (f(Q))- An equivalent statement is the 
intensities about Q = 0, 0, and Q = 0, 0, 2 are sym- 
metry related. This strong condition on the scattering 
pattern is in serious contradiction with the experimen- 
tally observed results. In contrast, if we consider an 
anisotropic Heisenberg pyrochlore model fEa. 1161 with fi- 
nite A) , we find that the lattice and spin degrees of free- 
dom do not force the scattering intensity to be identical 
about Q = 0,0,0 and Q = 0, 0, 2. For a model with full 
0(3) spin symmetry, the scattering profile is controlled 
by F"j_(q), Eq. El The significant difference between 
Eas. lT71 and H^l is the restoration of spin isotropy, i.e., 
the geometric factor defining the local (111) quantization 
axis, z a , is absent from Eg. 1171 Therefore, on purely sym- 
metry grounds, no (111) Ising model (i.e., Hamiltonian) 
with arbitrary distance dependent Jjj(rij) for Tb2Ti2 07 
will reproduce the experimental PM correlations shown 
in Fig. [3Ja) . Earlier works have recognized the need to 
consider more isotropic spin models for Tb2Ti207. In 
Ref. 0, the qualitative features of the PM scattering 
were reproduced from an isotropic structure factor for 
the nearest-neighbor Heisenberg pyrochlore AFM. Simi- 
lar results were obtained in Ref.|42by considering specific 
Q = spin structures on a cluster of two tetrahedra. 

To support the picture obtained on symmetry grounds, 
we have applied our MF formalism to two models for 
Tb2Ti20y: (1) a pyrochlore system with (111) Ising 
spins, and (2) a pyrochlore lattice with Heisenberg spins 
and finite (111) anisotropy. In both cases, the dipole- 
dipole interactions are evaluated with the Ewald method, 
see Appendix [0 For the (111) Ising description (i.e., 
A/| J\ > 1 and A/D > 1 and Eq. UHJ), we use Ji = 
-2.64 K and D = 0.48 K (Ref. 25). Our data are shown in 
Fig-IHb). Note that the scattering about Q = 0, 0, and 
Q = 0, 0, 2 is symmetry related, as predicted above, but 
is an intensity minimum. Monte Carlo data for this (111) 
Ising model agrees with our MF results, see Fig.|3J Monte 
Carlo simulations for Tb2Ti207 as a (111) Ising dipolar 
model were performed on a L = 4 lattice (N = 1024 
spins) at T = 5 K (T C MF ~ 1 K) with J x = -2.64K and 
D = 0.48 K, with the dipolar sum treated via the Ewald 
summation method, neutron-scattering data (as deter- 
mined by Eq. ^2 (Ref- 32)) were collected after 5xl0 7 
Monte Carlo steps per spin for both equilibration and 
measurement stages, and are shown in Fig. [3] The inten- 
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FIG. 2: (Color online) Paramagnetic scattering in the 
(hhl) plane for Tb2Ti2 07 : (a) Experimental paramagnetic 
scattering— (T = 9K), maximum intensity at Q = 0, 0, 2, (b) 
MF model of Tb2Ti2C>7 treated as a (111) Ising pyrochlore 
(T = 1.5T C MF ), no intensity at Q = 0,0,2, (c) MF model 
of Tb2Ti2C>7 treated as an anisotropic Heisenberg pyrochlore 
(T = 1.5T C MF , A = 20 K), maximum intensity at Q = 0, 0, 2. 



sity minimum at Q ~ 0, 0, and Q = 0, 0, 2 supports the 
above mean-field results and symmetry arguments. For a 
Heisenberg model with finite anisotropy (i.e., A/| J\ > 1 
and A/D > 1 and Eq. I16fl . our MF results are provided 
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in Fig |2I c) . With an anisotropy strength of A = 20K 
(i.e., A/D 41.7), we achieve good qualitative agree- 
ment with the experiment. The region around Q = 0, 0, 2 
has the strongest scattering with reduced intensity near 
Q = 2, 2, and the interconnecting regions. If we turn 
off the finite anisotropy (A — 0), i.e., an isotropic Heisen- 
berg model with long-range dipoles, the dominating scat- 
tering remains about points Q = 0, 0, 2 and Q — 2,2, 0, 
but there is increased intensity along the bridge regions in 
q-space connecting these points. Finally, in the absence 
of dipoles and A = one has the near neighbor AFM ex- 
change Heisenberg model, where the scattering intensity 
forms a network of interconnected triangles with equal in- 
tensity about Q = 0, 0, 2 and Q = 2, 2, fRefs. l4lll6l65h . 
Hence, a partial restoration of the spin isotropy is suffi- 
cient to place scattering about points Q — 0, 0, 2 and 
Q = 2, 2,0 in q-space, but to achieve good qualitative 
agreement with the experimental intensity profile dipolar 
interactions are necessary as is a finite single-ion contri- 
bution to the Hamiltonian. 
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FIG. 3: (Color online) Monte Carlo results for paramagnetic 
scattering in the (hhl) plane for Tb2Ti207 as a (111) Ising 
AFM, T = 5. OK, N = 1024 spins. Note that there is no 
intensity about Q — 0, 0, 2. Dipoles were treated via the 
Ewald method. 



IV. DISCUSSION 

A. Tb 2 Ti 2 7 

The prediction of an (111) Ising model for Tb2Ti2C>7 
is that of an AFM long-range ordered state in which all 
spins point either in or out of the unit tetrahedra at 
T w 10° K. However, this prediction is not realized ex- 
perimentally. Another problem, and possibly even more 
important, with a (111) Ising model for Tb2Ti207, is 
that the mean-field PM correlations, see Fig.^b), do not 
agree with the experimentally observed results, Fig.|2a). 



In Section IHII and in Appendix [D] we demonstrated on 
symmetry grounds alone that no pyrochlore model with 
(111) Ising moments could reproduce the basic features 
of the experimental PM scattering, i.e., strongest inten- 
sity centered about Q = 0,0, 2, lower intensity about 
Q = 2,2,0, and minimum intensity at the zone center 
Q = 0,0,0. These symmetry arguments were also sup- 
ported by MF and MC calculations of the elastic neutron- 
scattering cross section, Eq. ^] 

In Section II I II and in Appendix [DJ we were also able 
to demonstrate on symmetry grounds that a Heisenberg 
pyrochlore model of Tb2Ti207 would allow for elastic 
scattering about Q = 0, 0, 2 while at the same time per- 
mit no scattering about Q = 0,0,0. A MF calculation 
of the PM neutron-scattering in the (hhl) plane for an 
anisotropic Heisenberg model with long-range dipoles, 
Eq. ^ shows good agreement with the experimental scat- 
tering pattern (strongest intensity about Q = 0, 0, 2 with 
reduced scattering at Q — 2, 2,0 and along the bridges 
between these points, see Fig. EIc)). Both the finite 
anisotropy (A/|J| > 1 and A/D > 1) and long-range 
dipoles are necessary to achieve a quantitative match 
with the experiments. Reducing either the single-ion 
anisotropy to the isotropic limit (A = 0) or the range of 
the dipole-dipole interactions in the model reduces this 
agreement by altering the ratio of scattering intensity 
between Q = 0,0,2 and Q — 2,2,0. However, even a 
dramatically simplified model which would have dipoles 
cutoff at the first nearest-neighbor, r c = 1, would still im- 
prove the picture provided by a nearest-neighbor AFM 
exchange only Heisenberg model presented in Ref. ItH 
Therefore, and foremost in our argument, a restoration 
of the spin isotropy is absolutely necessary to place para- 
magnetic scattering about the q-space points Q = 0, 0, 2 
and Q = 2,2,0. The dipolar interactions are then im- 
portant for shifting (i.e., redistributing) the scattering 
intensity from Q — 2,2,0 to Q — 0,0,2. Intermedi- 
ate regions between these two points also experience a 
reduction in scattering. In terms of the underlying soft- 
mode spectrum, the q = eigenvalues and eigenvectors 
control the scattering at Q = 2, 2, to Q = 0, 0, 2. A 
shift in intensity from Q — 2,2,0 to Q = 0,0,2 signals 
a PM spin structure that prefers to lie in xy-plane (i.e., 
neutron-scattering at Q = 0, 0, 2 comes from spins with 
components perpendicular to this direction). 

Switching to a model with fully isotropic Heisenberg 
spins (as in Ref. l4l|) restores all the spin symmetry in 
the paramagnetic limit. This picture is dramatically in- 
consistent with the experimentally determined single-ion 
structure of Tb 3 + (J = 6, 7 F 6 ) in Tb 2 Ti 2 7 , where a 
ground state doublet is separated from the first excited 
doublet by an anisotropy gap of 20 K, close to the 9cw 
temperatur e) 24 ! 25 . Therefore, a restoration of the full spin 
symmetry at T < 20 K seems an unlikely explanation for 
the PM scattering at 9 K. 

The current MF approach does not allow the single- 
ion properties to be systematically considered, but a RPA 
calculation doesSi. By retaining only the simplest energy 
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level structure in the Tb 3+ wave function, the ground 
state doublet and the first excited state doublet, one can 
relax the strict (111) Ising constraint on the spins in a 
controlled approximation. Within the RPA, fluctuations 
out of the ground state and into the first excited state lev- 
els are equivalent to a fluctuating canting of spins away 
from the Ising geometry. In this case, the lowest or- 
der fluctuations from the strict (111) Ising limit yield 
qualitative agreement with experiment for the param- 
agnetic spin-spin correlations 66 . Others have also pro- 
posed a simple relaxation scheme of the strict (111) Ising 
directions 4 ^. Theoretical work remains to be done to ex- 
plain the failure of Tb2Ti2 07 to order at a temperature 
of 1 K, and why it remains paramagnetic down to 50 mK 
(Refs.Eland||3). 



B. General Discussion: Avenues for Other Studies 

We now briefly discuss some puzzling experimental re- 
sults for a few highly frustrated magnets. We note that 
the present mean-field formulation for the structure fac- 
tor, S(q), could provide valuable insight on the devel- 
opment of magnetic correlations out of the PM regime 
for each of these systems. The first, very paradoxi- 
cal, system is the antiferromagnetic GdaGasOi^ garnet 
(GGG). This material, where Gd 3+ is the magnetic ion 
with a spin S = 7/2, consists of two sublattices of in- 
tertwined spirals of corner-sharing triangles. For clas- 
sical Heisenberg spins coupled by nearest-neighbor an- 
tiferromagnetic exchange, each spiral on a garnet lat- 
tice structure should display a thermally-induced spin- 
nematic order-by-disorder transition according to work 
by Moessner and Chalked. Some precursors of spin 
coplanarity in GGG may have recently been observed 
in Mossbauer experiments^. In GGG, however, dipo- 
lar interactions are approximately 50% of the strength 
of the exchange interactions for nearest-neighbors and 
are, consequently, a sizable perturbation to contend with 
in this system 4 ^. In zero applied magnetic field, specific 
heat, magnetic susceptibility and nonlinear susceptibility 
measurements on GGG strongly suggest that this ma- 
terial undergoes a spin-glass transition around 140 mK 
(Refs. l45l and Hit) . However, the nonlinear susceptibil- 
ity, Xah measurements indicate that the spin-glass tran- 
sition in this material is unusual in that Xni exhibits 
two maxima 46 . In contrast to bulk measurements^^ 6 -, 
neutron-scattering experiments on powder samples of iso- 
topically enriched 160 Gd (natural Gd has a huge neu- 
tron absorption cross section) indicates the development 
of spin-spin correlations at approximatively 140 mK, ex- 
tending to a length scale of ~ 100A (Refs. E3 and Esfr . 
It is unclear at present whether or not the development 
of extended spin correlations in GGG at ~ 140 mK is 
an intrinsic effect or is due to material impurities and/or 
defects (e.g., Gd 3+ magnetic ions at Ga 3+ sites 67 ). An- 
other interesting system is the Gd2Ti2 07 pyrochlore an- 
tiferromagnet, where Gd 3+ is also the moment-carrying 



species. In Gd2Ti207, the dipolar interactions are ap- 
proximately 20% of the strength of the exchange in- 
teractions for nearest-neighbors and is here, just as in 
GGG above, an important perturbation 4 * 5 ^. Palmer 
and Chalker argue that the ground state consists of a 
fully ordered structure where each tetrahedral unit cell 
has an identical (zero total magnetic moment) spin con- 
figuration (a so-called q or d = structure) 5 -. Recent work 
has confirmed that this ground state is extremely robust 
against quantum fluctuations 53 . However, recent exper- 
iments on Gd2Ti207 are rather puzzling and appear in- 
consistent with Palmer and Chalker's work. Specifically, 
neutron diffraction measurements on 160 Gd isotopically 
enriched powders find a partially ordered phase with one 
disordered sublattice and propagation vector q = | , | , | 
at T = 50 mK (Ref . 6) , hence incompatible with the pre- 
dictions of Palmer and Chalkes 5 -. Specific heat measure- 
ments find strong evidence for two transitions at T = 0.7 
K and at T = 1.0 K (Ref. Recent mean-field cal- 

culations find evidence for a two-step magnetic ordering 
in this system. 68 69 So for Gd2Ti207, there also exists 
a complex behavior as signaled by thermodynamic mea- 
surements, theoretical predictions and neutron-scattering 
results. 51 Finally, the Yb2Ti207 pyrochlore ferromagnet 
is also puzzlin g 26 i 52 . There, neutron-scattering results, 
muon spin relaxation and Mossbauer experiments sug- 
gest a ground-state that lacks long-range magnetic order 
while there exist good evidence from the Mossbauer data 
that a first order spin freezing transition occurs around 
Tf ~ 0.24 K. Meanwhile, elastic neutron-scattering re- 
sults reveal the development of nontrivial spin-spin cor- 
relations as Tf is approached from above. In this sys- 
tem too, it is possible that long-range dipolar interac- 
tions may play some role due to the contribution of 
the I J, mj) — 1 7/2, ±7/2) eigenstate within the effective 
S e g — 1/2 ground state doublet. 26 



V. CONCLUSIONS 

In conclusion, we have demonstrated on symmetry 
grounds, through MF calculations, and MC simulations 
that the experimentally measured PM elastic neutron- 
scattering in Tb2Ti207 is inconsistent with a (111) Ising 
pyrochlore spin structure. From the qualitative agree- 
ment obtained using an anisotropic Heisenberg model, 
we argue in favor of a more isotropic effective spin model 
to describe the low energy phases of Tb 2 Ti 2 07. 

Finally, we have discussed the usefulness of a combined 
mean-field theory and Ewald method approach to study- 
ing geometrically frustrated magnets with long-range 
dipole-dipole interactions in the paramagnetic phase. 
This approach could be applied in general to any geo- 
metrically frustrated system. The zero field picture of 
Gd 3 Ga 5 0i2 (GGG) is particularly interesting because 
the low temperature phase remains an unresolved issue 
that entails unraveling the physics of competing exchange 
and long-range dipole-dipole interactions in a garnet lat- 
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tice environment 45 i 46 . 47 i 48 i 49 i 70 . 
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APPENDIX A: NEUTRON SCATTERING IN 
THE GAUSSIAN APPROXIMATION 



where J^{i,f) is defined by Eq. and a fictitious field 
term has been added. The field term, with \h%\ = \h\, is 
removed from the final equations by taking hf — > 0. We 
note that indices (a, b) = 1,2,3,4 label the sublattices 
and (u, v) = 1,2,3 label the spin components. 

In VMFT, an approximate free energy as a function of 
a trial density matrix, p, is formed, 

F p = Tr{pH H } + TTr{plnp} = (H H ) p +T(lnp) p , (A2) 

where Tr represents a trace over spin variables. F p is 
variational and defines an upper bound to the actual free 
energy, i.e., F p > F. The best functional from for trial 
density matrix is obtained by minimizing F p with respect 
to the parameters of p. For a system of N particles, the 
MF form of the A^-body density matrix is given by a 
product of single particle density matrices, 



p({sn) = l[Pi(S!) 



(A3) 



As we noted in the Introduction, other au- 
thors have discussed MFT and its application to 
magnetismSi 5 ^ 5 ^ 5 ^^ 7 . and the Ewald method for mag- 
netic dipole a 59 ! 60 ! 61 . Our purpose is to combine the tech- 
niques of MFT, developed here, with the and Ewald pro- 
cedures for magnetic dipoles in q-space, developed in Ap- 
pendix [CI so they can be readily applied to other prob- 
lems of highly frustrated rare-earth magnets. 

In this appendix we provide a detailed derivation of 
the mean-field equations for the elastic scattering cross 
section for pyrochlore spin systems. Our derivation is 
performed for Heisenberg spins, with a finite local (111) 
anisotropy, in order to broaden the appeal of the results. 
Connections to (111) Ising systems, infinite local (111) 
anisotropy, are noted at appropriate points. As men- 
tioned in Section [H] the MFT is developed via a varia- 
tional approach (VMFT), and in general this approach 
applies to a large array of statistical models with arbi- 
trarily complex order parameter In this work, the 
method reproduces the Gaussian approximation (GA) 
of the Landau free energy. The VMFT described here 
has been used by other o 39 i 55 i 56 i 57 , and we provide a de- 
tailed presentation here to clear up the notational incon- 
sistencies that appeared in some of these previous works. 
We also wish to provide for comparison with the RPA, 
which allows for a more controlled relaxation of the (111) 
Ising restrictio n 6 ?- 66 , that experimental evidence suggests 
is needed for Tb 2 Ti 2 0^i^. 

We begin with the model Hamiltonian of Eq. ^ Hh- 
The conventions for the spin vectors, Sf, unit vectors, 
n u and z a , and the description of the pyrochlore lattice 
in a rhombohedral basis with a four atom unit cell are as 
described in Section^ Therefore, our starting Hamilto- 
nian for MFT is given by, 



H H = - 



b,v 
3 



The single particle density matrix, p", is treated as a 
variational parameter that is subject to the constraints 



Tr{p?} = 1 , 
Tr{pfSf } = ml , 



(A4) 
(A5) 



which keep the internal energy constant, i.e., Tr{pH} = 
C. Here, is a vector order parameter, for (111) Ising 
spins one has the scalar equivalent, Tr{p" a? } — mf. In- 
corporating the constraints into the expression for F p 
gives, 

F p = Ti{pH H }+TTr{p\np} (A6) 

- ttt{J2$(p1-i)} 

i,a 

T Tr{^(p»5f - ml) ■ A°;} , 

where £f and Af are the Lagrange multipliers for the 
constraints of Eq. IA4I and Eq. IA5I respectively. In min- 
imizing F p with respect to pf, one has the following re- 
sults, 

-JLtt{ P H h } = , 

Spf 

-^Tr{plnp} = Tr{lntf} + Tr{l} , 
^Tr{££>-} = Tr{&}, 

i i,a 

^-T Tr{^^S? • A®} = Tr{Sf ■ A?} . 

" l i,a 

The optimum form for the density matrix is found by 
solving SFp/Sp? = 0, 



i,j a, b u,v 



(Al) 



(A7) 
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where Cf = e^" 1 = l/Tr{e A ? s ?} follows from Eg. 1X1 
Evaluating the trace in C°" for both Heisenberg and Ising 
spins we obtain 

(2tt) 3 / 2 

zf = hk^ 1112 ^ (A8) 



4tt 



and 



^sinh(|A?|) 



Z?=2cosh(A?) 



(A9) 



respectively, where i"i/2(|A"|) is a modified Bessel func- 
tion. The variational local density matrix is now written 
as 



(A10) 



and is used, along with the constraints of Eqs. IA4I and 
5J to rewrite the variational free energy, 



F ? = -^EE^^^' (A11) 

i,j a,b u,v 



We want an expression for the free energy to quadratic 
order in m°. This means that one must expand 
\nZ°"(Af) and then express A" as a function of mf. 
From the series representation of Z°"(Af) followed by the 
series expansion of ln(l — x), one has 



In Zf « InCi 



la\2 



2n 



(A12) 



where C\ is a model dependent constant and n = 1, 3 for 
Ising and Heisenberg spins, respectively. Using Eq. IA5I 
one obtains the expression 



mf = A 



al2\ 



l/2{ 

= i?<|coth(|A»|) 



1 



A? 



(A13) 
(A14) 



for Heisenberg and 

m? = tanh(A?) (A15) 
for Ising spins. To first order, we have 

m? = iAI , (A16) 

and 

m? = A? . (A17) 



Using Eas. lA12IA17l we can write the MF free energy to 
quadratic order in the order parameter, 

t,j a, 6 

- h T™T - TpN ceU In Ci , (A18) 

where p = 4 denotes the size of the basis (sublattice). 
As a side note, the Lagrange multiplier A® can be inter- 
preted as an effective mean-field interacting with a local 
moment. Minimizing F p: of Eq. IA11I with respect to 

a. ii. i 



T Q 



where 



T ' 



(A19) 



(A20) 



is the u-th component of the effective field at site (i, a). 
We next exploit the fact that the pyrochlore lattice has 
the underlying symmetry of a fee lattice by defining the 
Fourier transforms, 



E 



(A21) 



= T^-E^(9)e ig -< , (A22) 



where iV cc n is the number of fee Bravais lattice points, 



-*ab pfl jjb 

ij 

Equations IA21I and IA22I applied to F p yield 



J2f denotes the position of a spin, and Rfj = Rf — Rj, 



b,v 



r P (T) = ^EEE" 1 "^"'"-^^)) 

q a,b u,v 

- EEE rt "'"' ! " m V r P lnC i >( A23 ) 

q a,b u,v 

where T P (T) = F p (T)/N cc \\. A transformation to nor- 
mal modes is necessary to diagonalize J^{q)- This is 
accomplished by the use of Eq. or 

4 3 

m r = EE^wr- ( A24 ) 

a — 1 ^l— 1 

where the Greek indices (a, /i) label the normal modes. 
U (q) is the unitary matrix that diagonalizes J(q) in the 
sublattice space with eigenvalues A(q), 



U\q)J(q)U(q) = X(q) , 



(A25) 



where in component form U®'^(q) represents the (a, u) 
component of the eigenvector at q with eigenvalue 

A"(g). The amplitudes of the normal modes are denoted 
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by {4>q} — {J2 a n^q'^ 1 }- Therefore, the MF free energy 
to quadratic order in the normal modes variables reads, 



so one has, 



q,a,n 



(A31) 

where a constant term has been dropped. The integral 
above is recast as a general Gaussian, 



where 

a,u 

Note that for the Ising case, the indices representing 
the spin components [u, v) and the corresponding modes 
(/i, v) are dropped from Eq. IA26I 

The neutron-scattering scattering cross section for un- 
polarized neutrons in the dipole approximation is given 
by Eq. EU (Refs. E2 and El, or 



da(Q) _ C[f(Q)]' 



(19. 



cell 



EE^- S ix)e ,Q '< , (A27) 



i,j a, 6 



where Q is the momentum transfer, i.e., Q = G + q, G 
is a fee reciprocal lattice vector, q is a wave vector in the 
first Brillouin zone, and C is a constant. The correlation 
function is between spin components perpendicular to 
the vector Q, 



(Si ± ■ Sjj_) 



E 



(nl-nl)(SrS b i v ) 



(A28) 



Q U Q V \ 



where = h u — (n u ■ Q)Q is strictly a geometric factor 
and Q = Q/\Q\. For Ising spins one replaces n u with z a 
from Table H] and Sf' u with erf. In order to proceed, the 
correlation function between spin variables, (Sf ± ■ Sj ± ), 
must be transformed to q-space by use of Eq. IA21I and 
then to normal mode variables by application of Ea. lA24l 
one arrives at 



EEEW^) (A29) 

q,q' a, 

U^{q)Uti{q')o- lq - R '- lqhR " • 



The correlation function {4>%^4>qi V ) can be calculated 
from a partition function defined in terms of the normal 
mode amplitudes. The general definition reads, 



Z = Tr{e 



-F P (T)/T 



} 



(A30) 



where T P {T) is given by Eq. IA26I and the trace is over 
all values of the normal mode amplitudes, 



/oo 
n > 
-oo „ „, ,. 



l -A4> 2 +Bcl> _ 



I2n b 2 



(A32) 



where A = (n ^) and B = h%». Therefore, the 

final form of the partition function is, 



z = n zc ^(q) 



(A33) 



n 



2tt 



1 (1/2) 



T -I 



J^' M | 2 /2("-A°(q)/T) 



The correlation function is now determined from deriva- 
tives of Z with respect to the fields h^v, 



yxvjfi^ - I d H_ 



(A34) 



Back substitution of the result from Eq. IA34I into 
Eqs. lA"29l lA"28l and then into Eq. \M7\ and finally im- 
posing the properties of the Kronecker delta functions 
leaves 



1 da{Q) 



A^cdi dn WWJ (n-A«(g)/T) 

ct,fj, a,b U,v V (A / 



x ^(g)^(-g)e lG " 

where we have used the identity 

^ e *(Q-q)-K lJe *(Q-<j)-r« 6 = ^^G-r" 6 _ 

We find it convenient to define the function, 
KJl) = ^^K^Y^ (A36) 

a , u 
a 

where £/£' a (q) = E^^m/uO?) ( and therefore F° x (g)) 
is a three component vector. The scattering cross section 
is written compactly as 



(A35) 



JV ceU dO 



a . / j. 
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which is Eq. ^] with n = 3. When considering (111) 
Ising spins, the arguments that follow Ea. lA29| still hold, 
but the indices for the spin components (u, v) and cor- 
responding normal modes (/x, v) are dropped from the 
equations. The final expression for the scattering cross 
section reads, 



in Eas. lA36l and lA39l is replaced by exp (iQ ■ r a ). Recall 
that momentum transfer and reciprocal lattice vector are 
related according to Q = G + q. The presentation found 
in Refs. y!J and |56l (as noticed recently by Kadowaki et 
alJ2i) unintentionally mixes these two conventions. 



N, 



1 da(Q) 



cell 



= C[f(Q)} 2 Y, 



(1-A«(g)/T) 



(A38) 



where F"(q) is a three component vector given by, 



[zlU a ' a (q)e 



iG-r a 



(A39) 



which is just Eq. [18J We note that Eqs. IA37I and IA38I 

are only valid for T > A^(g rd)/»i, where A^(q rd) = 
maxq{A max (<7)} is the critical eigenvalue, a global max- 
imum, which sets and defines the paramagnetic 
regime. We also point out that the lattice structure and 
spin symmetry are contained in the respective F±(q)- 
functions, Eqs. IA36I and IA39I these properties will be 
exploited in Appendix iDl when we discuss symmetry al- 
lowed scattering patterns for Heisenberg and Ising spins. 



Comment on the convention for the Fourier 
transform 



APPENDIX B: HIGH TEMPERATURE 
EXPANSION OF THE GAUSSIAN 
APPROXIMATION 

We demonstrate that the equations for the neutron- 
scattering cross section obtained from MFT at the Gaus- 
sian level, AppendixE] can also be formulated via a high 
temperature series expansion (HTSE) to lowest order in 
(3 = 1/T, where T is temperature in units of fee- In con- 
trast to VMFT, in a HTSE there is no appeal to any 
simplifying approximation that changes the character of 
the density matrix and imposes constraints that keep the 
internal energy, Tr{pH}, fixed. We follow our established 
convention of treating the general case of anisotropic 
Heisenberg spins while pointing out the specific differ- 
ences for (111) Ising spins when needed. The starting 
point is the Heisenberg Hamiltonian of Eq. 



(Bl) 



i,j a,b u,v 



Another convention for defining Fourier modes uses the 
Bravais lattice points, i?,, in the definitionsii, 



E 



m a,u e - lq - Rl 



and 



^uv(h j) 



A^ccll 



E 
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«9)e 



iq-Ri 



(A40) 



(A41) 



Hence, this convention differs from the one we employ 
by a simple phase factor, exp (iq ■ r ab ) for J£%{i,j)- The 
two approaches are equally valid and yield the same re- 
sults; however, there are a couple important differences 
in the above results when this alternate convention is 
used. First, the interaction matrix defined by the inverse 
of Eq. IA41I has complex entries. For nearest neighbor 
interactions between sites a — 1 and b — 2, we have 



<jl,2 



(i 



-iq-Ri 



(A42) 



where the factor of 1 arises because the sites a = 1 
and 6=2 are in the same tetrahedral basis unit, i.e., 
Ri,j = 0, but its symmetric equivalent has Ri.j ^= 0. 
Next, the definition for the scattering cross section, 
Eq. IA27I holds, but the expression for the correlation 
function between normal mode variables, Eq. IA29I con- 
tains the phase factor exp (— iq ■ Ri — iq' ■ Rj). Carrying 
the necessary steps through to an expression for scatter- 
ing cross section, one finds that the factor exp (iG ■ r a ) 



where J£%{i,j) contains both spin and coordinate de- 
grees of freedom and is defined by Eq. 

In the formula for the scattering cross section, Ea llll 
one must calculate the perpendicular correlation func- 
tion, 

( S i± ■ S k) = E("" • bUWSp) , (B2) 

u,v 

which is just Eq. El We express the correlation function 
as a series expansion in cumulantsS, 

00 (— R\ m 

{S a,u s b. V) = J2 1 — L_| — (Si' U Sj' V Hff) c , (B3) 



where (...) represents a thermal average with respect to 
Hh and (...) c represents the cumulant expansion of the 
spin operators and Hh- Cumulants are evaluated as a 
trace over the T = states, (...)o = Tr{...}/Tr{l}, where 
Tr{l} = N is a normalization factor, (4^)^ for Heisen- 
berg spins and (2) N for Ising spins, N is the total number 
of sites in the lattice. Therefore, the correlation function 
to lowest order in (3 is 



{S«' U S°' V ) « (S? U S$' V ) C - (3{S*' U S°' V H H ) C + ■■■ (B4) 

For both Heisenberg and Ising Hamiltonians, any non- 
zero contribution to the correlation function must have 
an even number of spin components per site (i.e., 
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(S^ ,u ) 2 ). The zeroth order term in j3 is obtained triv- 
ially, 

{S a,u s b,v )c = {S a,u s b,v )o = ^ / ^.^b gu,v ^ 

The first order contribution has two terms in the cumu- 
lant, 

{S a,u S b,v H)c = {s a,u s b,v H)o _ ^gb.v^^ ^ 

but the second does not contribute because (H) a oc 



where A"(— q) = A"(q). Using the results from Eas. IB8 



and IB9I the expression for the scattering cross section 
becomes, 

^^-(c/»)l«<M*i:£2>l..l) 



a,fi a,b u,v 



„iGr" 



(BIO) 



(Sf' u S^ l v ) = 0, spins in the Hamiltonian can not be 
at the same site. The first term yields the result 

(B5) 



In the high temperature limit (/3 — > 0), one can write 

-l 



1 + 



1 



l,m ^ £ 
2\ rrabi 



Therefore, Eg. IB4l can be rewritten as 



(SrS°> v ) « (l/n)S id 6 a ' b S^ + (P/n 2 )j:Z(i,j) (B6) 

We obtain an expression for the scattering cross section 
by substituting the result of Ea. lB6l into Ea. lB2l and then 
that result into Eq. ^1 we use the identity Q = q + G 
and the definition of the Fourier transform of J^{i,j), 
Eq. to arrive at, 



1 da{Q) 
N ce n dQ 



= (C/n){f(Q)} 2 Y,Y,^-^ ( B? ) 



a,b u,v 



x [ 5 a - b S u - v + -J^{-q) I e 
n 



iGr a 



where C is a constant. Note, a factor of l/iV co n has been 
absorbed into the expression of the Fourier transform 
of the correlation function, Eq. IB6I in order to remain 
consistent with our use of Eq. IA23I in the derivation of 
the mean- field neutron-scattering cross section, Eq. IA37I 
The interaction matrix, J^{— q), in Eq. IB 71 is not diag- 
onal. In order to calculate the differential cross section 
for the pyrochlore lattice, or any lattice with a basis, we 
must diagonalize J^{—q)- This is done with unitary 
matrix U(q), or U£'"(q) in component form. The first 
term on the right hand side of Eq. IB 71 follows directly 
from the definition of U(q), 

E K:»{-<l)U^{q) =1 = 5 a > b 5^ , (B8) 

a, fx 

where I is the appropriate identity matrix. The trans- 
formation of the second term uses Eq. IA25I and solves 
for J(q), 

J(q) = U(q)X(q)U\q) , 
or in component form, 



and the mean-field result is recovered, 



1 MQ) -cmQ)YY. lF ' Aq) 



N, 



n-/3A«(g)) 



(Bll) 



where F^ ± (q) is given by Eq. IA36I with U^(—q) = 
U^u(q)- For a (111) Ising model, the indices (u, v) and 
(/i, u) are dropped from our presentation, n = 1, and 
F£(q) is given by Ea. lAl9l 



APPENDIX C: EWALD EQUATIONS 

Here we treat the dipole-dipole term in J{q) via the 
Ewald method^. In MFT one works in the thermody- 
namic limit, so one has an infinite lattice sum. Within 
the Ewald approach, one recasts this infinite and condi- 
tionally convergent series as two finite absolutely conver- 
gent (rapidly converging) sums£2*£^. The application of 
Ewald's ideas to the Fourier transformed dipole-dipole in- 
teraction is equivalent to the method of long wave lengths 
presented in Ref. |HE 

The general expression of the Fourier transformed 
dipole-dipole lattice sum is 



A(q) 



E 



'E<< 



zq-R% 



(CI) 



where 



Aab _ 



3(n" ■ R^)(n v ■ i?° 



(C2) 



The conventions for indices and vectors are described in 
Section m and A(q) is a 12 x 12 matrix. The sum ^ 
is over all Rf b except the terms Rfj = 0. To implement 
the Ewald method, we rewrite Eq. IC2I in the following 
forma, 



(B9) 



>C = -(n t, -V s )(n"-V !e ) 



(C3) 
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Equation IC3I is used in Eq. IC1I and in terms of compo- 
nents one has, 



x=0 

(C4) 



AtM = -(n u ■ v,)(<i« • v.) ( E'SS: 



The goal of the Ewald method is to rewrite Eq. IC4I a 
conditionally convergent series, as two absolutely conver- 
gent series, one in real space and the other in reciprocal 
space (fc-space). We begin by writing the sum inside the 
brackets as a sum over all fl?Jr, the result is, 



-iq-Rf 



R?° - x\ 



6 ab (h u • V x )(n" .V,) n 
\x\ 



x=0 

(C5) 



Next, the definition of a Gaussian integral (also a gamma 
function identity^) , 



1 



-t R' 



dt 



is used rewrite the point source term, l/|fl" h — x\, in 



Eq. IC5I The Fourier transformed dipole-dipole lattice 
sum now reads, 



AtM = ~{n u ■ V x )(n v ■ V x ) / dt -=e-«- 
x j^e- t2 l R ? 



+ S ab (n u -V x )(n v ■ V x ) 



x 



(C6) 



x=0 



The integral in Eq. ICJ6I is divided into two regions, [0, a] 
and [a, oo). It is from this decomposition that the real 
space ([a, oo)) and fc-space ([0, a]) series will arise. We 
note that the series resulting from the [a, oo) integral will 
have a divergence at fly = 0, hence this term is treated 
separately. The range of integration is controlled by a; 
it has units of inverse distance and will play the role of 
a convergence parameter in the final series. Eauation lC6l 
now reads, 



where 



W^(q) = -(n"-V x )(n"-V x ) / dt^e-^ 

JO V 71 " 



fj2e- t2 \ R *h*\ 2 -*<i<<-*A ,(C8) 



K°M = -(n u • v*)(ft B • v.) / dt -r<r %q - x 



2|D ° 6 -x| 2 -j<j-(H° b -a 



,(C9) 



Y™ = 5 a \n u - V x )(n v -V x ) 

J__ _L r e -* 2 N 2 dA . (CIO) 

M J a ) x=0 

We treat the expressions for W^(q), X%*(q), and in 
succession. 

For W^(q), the sum inside the brackets is a periodic 
function in x. Therefore, it can be expressed as a Fourier 
series, 



Solving for one has 



9k=G 



An e 



-iGr a 



v \q-G\ 3 



(Gil) 



(C12) 



where G is a reciprocal lattice vector, i> is the volume of 
the unit cell, 



F{z) = / y sin(y) e z v dy 



4z 3 



-l/4 2 2 



(C13) 



and z = t/\q — G\ (Ref. 77). One now has the following 
identity for f(x), 



f(x) = £e ' " 



: <"$-x\ 2 -iq-(R$-x) 



4tT „ e -!G (r oi -o;) 
G 

Substituting Eq. IC14I into Eq. IC8I differentiating, and 
imposing the limit on x yields, 



2 = /"dt F(i/|g-G|) . 

I" JO 



(C15) 



The integral over [0, a] is readily performed by using the 
result from Eq. IG13I Therefore, the reciprocal space sum 
in the Ewald decomposition reads, 

W»(a\ = ^ v [n".(g-G)PMg-G)] 

wWJ « ^ |q-G| 2 



(C7) 



x e-l^ G l 2 / 4Q2 e- lG "-° 



(C16) 
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where the sum is over all reciprocal lattice vectors G. We 
note, however, the series for W£ b (q) has a nonanalytic 
term at G — when at the zone center, q = 0. This 
point is discussed below. 

The expression for X^ b (q), Ea. lC9l can be rearranged 
to obtain an identifiable integral. By reversing the sum 
and integral in Eo. IC9l we obtain, 

Kt(q) = -(n u ■ v x )(fi v ■ v x ) ' e ~^ 



x — — / dt e 
V 71 " J a 



lx=0 



(C17) 



The integral in Eq. IC17I can be expressed as a comple- 
mentary error function^, 

2 /*°° 2 

crfc(z) = — = / e~ x dx . 

A final form for X^(q) is obtained by first applying the 
differential operators in Eq. IC17I followed by taking the 
limit x — > 0, and then integrating to get, 



where 



(C18) 



CI ab I r>ab\ 



(C19) 



2a e~ a2 \ R "\ 2 erfc(a|J2#|) 



l^l 3 



ocyab f ryab 



4a 3 



(C20) 



6a 



*3 

3 erfc(a|-R?/ 

\R?j b \ 5 



Eauations lC18l(^20l form the real space sum in the Ewald 
decomposition of the dipole-dipole interaction. Note that 
the sum in Ea. lC18l is over all Bravais lattice displacement 
vectors Rij , with j fixed, except R^ = 0. Hence, the real 
space Ewald series is analytic everywhere. 

In treating the singular terms in Eq. IC10I one applies 
differential operators first to get, 



x, — >0 



(C21) 



(n u -h v ) 3(n u ■ x)(n" ■ x) 



+ (Slt(x) -S2%(x))} , 

where Sl a u b v (x) and S2 a u b v (x) are given by Ens. RTml and 
20| respectively, with Rf b replaced by a;. To evaluate 



the limit in Eq. IU2T1 one expands the exponential func- 
tion to 0(x 2 ) and the complementary error function to 
order 0(x 3 ). The result is the constant, 



yah = 4a 3 

30? 



[n ■ n )o 



(C22) 



Collecting the results of Eqs. l€T6l IC18IC20I and IU221 

we write the Ewald representation of the (/-dependent 
dipole-dipole interaction as, 



4a J 



(h u ■ n v )5 a ' h 



(C23) 



v ^ 

+ E' (S1%(R$) - S2£(Kf»)) e-^-< , 



where 



K uv {q-G) = 



(g-G))[n"-(q-G)} 

la-GI 2 



(C24) 



In our derivation of the Ewald equations there is no 
reference to a specific lattice structure. Therefore, the 
Ewald results encapsulated in Ea. lC23l hold for any lattice 
described by a set of translation vectors {Rf b }. Through 
the unit vectors n u (where local quantization axes can 
be treated by including a sublattice index, i.e., h a,u ), 
the freedom to define the spin symmetry (e.g., Heisen- 
berg, XY, Ising) has been ensured, too. For the work 
discussed in this article, we consider both Heisenberg and 
(111) Ising spins on the pyrochlore lattice. For Heisen- 
berg spins, A^ v {q) is calculated for all sublattices (a, 6) 
and spin components (u, v), the resulting A(q) is a 12 x 12 
symmetric matrix contribution to For (111) Ising 

spins, the sums over spin components are dropped and 
the local quantization vectors are substituted, z a . One 
calculates A ab (q) for all sublattices (a, b), resulting in 
A(q) a symmetric 4x4 contribution to J{q). For each 
pyrochlore model, A(q) is determined at every q-point 
in a mesh that covers the first Brillouin zone in the (hhl) 
plane. These matrices are stored and then used in the for- 
mation of J(q) to calculate the neutron-scattering cross 
section, Eq.^|or Eq.^J for a specified set of interaction 
parameters (i.e., J, D, A, T). Because A(q) is calcu- 
lated only for q in the first zone, the term K. uv (q — G) 
in Eq. IC23I is ill defined at q = G — 0. We discuss the 
small q behavior of the Ewald equations below. 

The parameter a used to divide the integral in Eq. IC6I 
functions as a convergence parameter in the Ewald sums, 
Eq. IC23I Although, the result of A ab v (q) is independent 
of the value of a, in practice one chooses an a so that both 
real and reciprocal sums converge rapidly. Note that the 
convergence of the real space sum, Eqs. IC18IC20I is en- 
hanced by a large value for a, while the convergence of 
the reciprocal sum, Eq. IC16I is improved for a small a. 
In choosing a convergence parameter, we followed Ref.l59l 
and set a = \Jir/v, where v is the volume of the unit cell. 
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For a pyrochlore lattice defined in the rhombohedral ba- 
sis with a cubic cell size of a, we used v = a 3 /A. The real 
and reciprocal space sums converged at about the same 
rate for this value of a. We obtained similar results for 
A^ v (q) using a — ^tt/2v and a = y/2w/v. Our Ewald 
results were checked by comparing the maximum eigen- 
values of A(q) to those generated from a direct lattice 
sum of A(q) out to some cutoff distance r c . Comparisons 
were done for the bcc and fee lattices. We also performed 
tests of our Ewald equations for the pyrochlore lattice by 
calculating the soft-mode spectrum of A(q) in the spin- 
ice regime, e.g., D = 1. Ewald results along (001) in the 
first Brillouin zone were compared to calculations with 
the dipolar sum cutoff at different maximum separation 
distances r c . The cutoff results approach the Ewald re- 
sults as r c increases. This spectrum of eigenvalues agrees 
well with the spectrum generated from a direct lattice 
sum for A(q) with a cutoff distance of r c = 1000, Fig. 6 
in Ref. |7£| The Ewald method eliminates the ripples in 
the soft-mode spectrum of A(q) by effectively taking the 
range of interaction to infinity. 



The reciprocal space sum in Ea. lC23l has a nonanalytic 
term at the point q = in the first Brillouin zone. If wc 
consider the G = contribution to Eq. IC16I we have, 



W£{q 



(C25) 



In the limit of small q the exponential is expanded to 
yield. 



W£{q,G = Q) 



4tt (n« • q)(h v ■ q) A |g| 



where in the limit q — > the value of 



4tt (h u ■ q)(n v ■ q) 



Aa 2 ) ' 
(C26) 



depends on the direction in which one approaches the 
zone center. The nonanalytic term can be related to the 
macroscopic field of the dipoles and is shape dependent, 
see Section 30 in Ref. |5^ We drop this term from our 
calculation to obtain a completely smooth spectrum all 
the way to q — 0, 0, 0. The physics of spin-ice is not af- 
fected by this omission because all modes contribute to 
the PM scattering with q = 0, 0, 1 going critical at T C MF . 
The case of Tb2Ti2C>7 is more subtle because it is the 
q = 0,0,0 soft mode that goes critical. However, our 
focus here is not the ordered state of Tt)2Ti207, where 
a q = ordered state is expected for a pyrochlore AFM 
with either (111) IsingSiSSiS or HeisenbergSi spins. In- 
stead, we are concerned with understanding the physics 
in the paramagnetic regime of this system as a first step 
toward unraveling the mystery surrounding the failure of 
Tb 2 Ti 2 7 to order at 50 mK. 



APPENDIX D: SYMMETRY EXCLUDED 
SCATTERING 

The paramagnetic neutron-scattering spectrum of 
Tb2Ti207 in the (hhl) plane contains a strong but broad 
region of intensity about Q = 0, 0, 2 with no discernible 
correlations near the zone center, Q = 0,0,0 (Ref. ^J). 
In this appendix, we put forward arguments based only 
on the structure of the lattice and the symmetry of spin 
space to demonstrate that the PM scattering intensity 
profile described above can not be realized by (111) Ising 
spins on the pyrochlore lattice but is allowed if the spins 
are Heisenberg like. 

For the Ising pyrochlores, the map of scattering inten- 
sity is determined by the function F"(q), Eq. ^3 which 
contains only information on the symmetry of the lattice 
through the eigenvectors, U a - a (q), and the phase factor, 
exp (tG ■ r a ), and the symmetry of spin space through 
the local quantization axis, z a . We consider a unit tetra- 
hedron with scattering vectors Q restricted to the (001) 
direction. To handle the situation near the origin, we 
express all Q as a small displacement from a reciprocal 
lattice vector, i.e., Q = G + q = 0, 0,£ + 0, 0, 5, where 
< S < 1, I is an integer, and a factor of 2tt/o, is im- 
plied. The term 0, 0, S falls in the first zone and, there- 
fore, determines the eigenvalues and eigenvectors. Using 
the values for r a and z a defined in Table |U we write, 

F2(o,o,e+S) = (u 1 ' a (s)-u 2 ' a (s)) (di) 



V3 

(1,-1,0) 
V3 



(U 4 > a (6)-U 3 ' a (6))~ l 



e" 2 . 



Note that the projections of the spins onto the plane 
perpendicular to the direction of Q sum to zero, (i.e., 

zy + Zj 2 ^ + z± + = 0). For wave vectors Q = 0,0,6 
and Q = 0, 0, 2 + 6 one has the following, 

(M,o) lTTl<a 



F±(P,0,6) 



V3 



(U 1 ' a (S)-U 2 ' a (S)) (D2) 



(!> ^°\u 4 ' a (S)~U 3 ' a (S)), 



and 



F1(0,0,2 + S) 



V3 

(1,1,0) fTTl 



V3 



U lta {S) - U 2 ' a (6)) (D3) 



(!, *'°\ u 4 ' a (6)-U 3 ' a (6)). 
v3 



The modulus squared of these two functions, e.g., the 
numerator of the scattering cross section, yields the same 
numerical result, 

|F«(0,0,,S)| 2 H*T(0,0,2 + <5)| 2 

= \{{U ha {S) - U 2 > a (5) - U 3 ^ a (S) + U^ a (5)) 2 

+ (U l ' a (5) - U 2 ' a (S) + U 3 > a (5) - J7 4 '"(<5)) 2 } . 
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This means the scattering cross section, given by Eq. 
in the limit 5 — ^ 0, is the same (or exactly correlated) 
for Q = 0, 0, and Q = 0, 0, 2, absent the magnetic form 
factor (f(Q)). Therefore, the paramagnetic scattering of 
Tb2Ti2C>7 can not be generated by a model with Ising 
spins (infinite local (111) anisotropy). 

In the case of Heisenberg spins with finite single-ion 
anisotropy, we consider the function F" ± (q) is given by 
Ea. 1171 Again, restricting ourselves to wave vectors along 
the (00/) direction, we have the general result 

F? t± (0, 0,£ + S) = (S) + (6) (D4) 
+ «'!W+<l(5))e^ , 

where E/£"(<f) = {U U«;™(5), 0). For Q near 0, 0, 
and 0, 0, 2, we obtain the following two forms, 

*£l(0,0,*) = U^(S)+Ul : a ± (S) (D5) 

+ 00 + 

and 

F« ± (0,0,2 + £) = UfcW + UfciS) (D6) 



Taking the modulus squared we get, 

|F" (0,0,<5)| 2 = (A 2 + B 2 ) + (C 2 + D 2 ) + 2(AC + BD) , 

(D7) 

and 

|F^'"(0,0,2+,5)| 2 = (A 2 +B 2 ) + (C 2 +D 2 )-2(AC+BD) , 

(D8) 

where 

A = UtZ(6) + U*<«(5) , 

B = U^(S) + U^(5) , (D9) 

C = U^(5) + U^(8) , 

D = U*$(6) + U*$(S) . 

Equations ID7I and ID8I are not strictly equivalent. Hence, 
it is possible to have paramagnetic spin-spin correlations 
about Q = 0, 0, 2 while intensity about Q = 0, 0, is 
suppressed. This result puts on a firm theoretical footing 
the need to describe Tb2Ti2C>7 by a three-component 
Heisenberg model with finite anisotropy. 
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